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Employing efficient diagonalization techniques, we perform a detailed quantitative study of the 
regular and chaotic regions in phase space in the simplest non-integrable atom-field system, the Dicke 
model. A close correlation between the classical Lyapunov exponents and the quantum Participation 
Ratio of coherent states on the eigenenergy basis is exhibited for different points in the phase space. 

It is also shown that the Participation Ratio scales linearly with the number of atoms in chaotic 
regions, and with its square root in the regular ones. 

PACS numbers: 


I. INTRODUCTION 

The non-equilibrium dynamics of isolated quantum 
many-body systems is a fundamental problem where rel¬ 
evant progress has been achieved, and many challeng¬ 
ing questions related with thermalization remain open 
(see mia and references therein). Altland and Haake 
mill have demonstrated that effective equilibration can 
occur for unitary dynamics under conditions of classical 
chaos, showing that the evolution equation of the Husimi 
function is of the Fokker-Planck type in the Dicke model 
[5]. The transition from ergodic to non-ergodic behav¬ 
ior in integrable many-body systems with a weak non- 
integrable perturbation has been quantified employing 
the average over an energy shell of the inverse participa¬ 
tion ratio between the eigenstates of the integrable and 
the total Hamiltonian. The delocalization is associated 
with the thermal behavior of the system[6l [7], however, 
in the case of the thermodynamic limit a further analysis 
is required [T]. 

The presence of chaos and its relation with the 
quantum phase transition (QPT) in the Dicke model 
were studied by Emary and Brandes [8] employing 
the semiclassical Poincare sections, the nearest neigh¬ 
bor distributions of eigenenergies and their similitudes 
with the Wigner or Poisson distributions. Employing 
the Holstein-Primakoff mapping, delocalization in phase 
space and macroscopic coherence are suggested as gen¬ 
eral features of the onset of chaos. Employing the efficient 
coherent basis (ECB) it has recently been possible to ob¬ 
tain exact numerical solutions both for the ground state 
and for a relevant part of the energy spectrum, 
both in the normal and superradiant phases, for a large 
number of atoms A/". With these tools we were able to ex¬ 
tend the above mentioned study to different energies and 
coupling constants. Classical and quantum qualitative 
signals of regularity and chaos were analyzed using, re¬ 
spectively, Poincare sections and Peres lattices, as well as 
the Anderson-Darling parameter as a quantitative tool to 
identify chaos in the quantum energy spectrum HIHISl, 
extending the pioneer work of Emary et al. [8]. 

In this work we analyze several regions in the phase 


space of the Dicke model, for different energies and an 
atom-field coupling strength, where the classical chaos 
conditions required in laia, characterized by their Lya¬ 
punov exponents, effectively occur, and exhibit that they 
are intimately associated with the Participation Ratio of 
the corresponding coherent state on the eigenenergy ba¬ 
sis. Following a seminal idea by Haake [14], it is shown 
that a large Participation Ratio of the coherent state, 
which quantifies its delocalization, provides a quantum 
measure of chaos in each point of the associated classi¬ 
cal phase space. The Participation Ratio of the coherent 
state is closely related with the Husimi function. As a 
by-product, we show that the Husimi functions for dif¬ 
ferent eigenenergy states in an energy window reproduce 
regular or chaotic Poincare sections, providing graphical 
and qualitative support for the previous findings. 

The article is organized as follows. In section 2 we 
describe the Dicke Hamiltonian and its classical limit. 
In section 3 we present some basic results characteriz¬ 
ing the classical and quantum dynamics of the Dicke 
model: Poincare sections and Lyapunov exponents, and 
the Husimi function of individual energy eigenstates. In 
section 4 we employ the Participation Ratio of coherent 
states in order to quantify chaos in phase space using the 
eigenstates of the Dicke Hamiltonian. Finally, we expose 
our conclusions. 


II. DICKE HAMILTONIAN AND ITS 
CLASSICAL LIMIT 

A. The Dicke Hamiltonian 

The Dicke model describes the simplest non-integrable 
atom-field system, exhibiting quantum chaos [HI nHi nnj. 
The Hamiltonian has three terms: one associated to the 
monochromatic quantized radiation field, a second one 
to the atomic sector, and a last one which describes the 
interaction between them. With h = I, it reads 

H]j = ujo^ a + ujqJzF 
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The parameter 6 allows to switch between the tradi¬ 
tional Dicke model (6 = 1) and its integrable approxima¬ 
tion (via the rotating wave approximation), the Tavis- 
Cummings model {S = 0). Here, the frequency of the 
radiation mode is cj, associated with the number opera¬ 
tor a^a. For the atomic part cjq is the excitation energy, 
while J^, J+, J_ are collective atomic pseudo-spin op¬ 
erators which obey the SU( 2 ) algebra. It holds that if 
j{j + 1 ) is the eigenvalue of = J‘^ ^ Jy ^ then 
j = M/2 (the pseudo-spin lenght) defines the symmetric 
atomic subspace which includes the ground state. The in¬ 
teraction parameter 7 depends principally on the atomic 
dipolar moment. Besides, commutes with the parity 
operator H, 


n = with A = a^a + -h j. ( 2 ) 

The eigenvalues A = n + m + j of the A operator are the 
total number of excitations, where n is the number of 
photons and riexc = m j the number of excited atoms. 
When (5 = 0, the Hamiltonian commutes with A, hence 
in this case it is integrable. For the rest of this arti¬ 
cle we limit ourselves to the ^ = 1 case which includes 
the anti-resonant terms, the non-integrable case or Dicke 
model proper. As mentioned before, one of the most rep¬ 
resentative traits of the Hamiltonian is its second-order 
quantum phase transition (QPT) in the thermodynamic 
limit 113120], a paradigmatic example of quantum col¬ 
lective behavior [21]. When the atom-field interaction 
reaches the critical value 7 c = + (5), its ground 

state goes from a normal (7 < 7 c), with no photons and 
no excited atoms, to a superradiant state (7 > 7 c), where 
the number of photons and excited atoms becomes com¬ 
parable to the total number of atoms in the system, i.e. 
a macroscopic population of the upper atomic level. 

Despite its simplicity the Dicke Hamiltonian remains 
as a model of great theoretical and experimental interest. 
A mean-field description of the ground state allows to ex¬ 
tract the critical exponents for the ground state energy 
per particle, the fraction of excited atoms, the number of 
photons per atom, their fluctuations and the concurrence 
[HIIISIES], however, around the QPT it has a singular 
behavior [23ti25] . Analytical expressions for its eigenen- 
ergies have been reported [26H3Q] . Another important 
feature in the Hamiltonian is the excited-state quantum 
phase transitions (ESQPT) [31]|32], manifested as a sin¬ 
gularity in the level density, order parameters, and wave 
function properties [33l [34] . They could have important 
effects in decoherence m and in the temporal evolution 
for quantum quenches [36]. It is strongly suggested that 
the relation between the ESQPT and chaos is dependent 
on the system m- 

While the Dicke Hamiltonian was designed to describe 
a system of Af two-level atoms interacting with a single 
monochromatic electromagnetic radiation mode within a 
cavity [5], it can also be employed to describe a set of 
JV qubits from quantum dots, Bose-Einstein condensates 
or QED circuits [38ll4T| . interacting with a bosonic field. 
It is worth to mention that the experimental observation 


of the super radiant QPT in a BEC system described 
by a Dicke-like Hamiltonian [42] has attracted renewed 
interest in its study. 


B. The classical Hamiltonian 

An effective classical Hamiltonian can be obtained 
employing Glauber and Bloch coherent states for the 
bosonic and pseudo-spin sector, respectively. In the case 
of the Dicke model this is a natural choice given the alge¬ 
braic structure of the degrees of freedom. Its dynamical 
properties can be described by the temporal evolution of 
this coherent state product, assuming the system remains 
in it [43] . The Glauber and Bloch coherent states for the 
bosonic and pseudo-spin sector, respectively, are defined 
as 


la) = 



(3) 


In order to obtain the effective classical Hamiltonian, we 
calculate the expectation value of the Hamiltonian opera¬ 
tor in the coherent state product [44] . As the dynamical 
description requires canonical variables, they are built 
from the complex parameters 2 : and a. Eor the Glauber 

parameter we consider a = + w) with q and p 

real values, whereas for the Bloch parameter, the stere¬ 
ographic projection ^ = tan(^/ 2 )e*'^ provides a set of 
canonical variable by considering jz = {jz/j) = — cos^ 
and = arctan(j^/jia,), where 0 and 0 are spherical angu¬ 
lar variables of a classical vector j = {jx^jy^jz) {\j\ = j) 
with 0 measured respect to the negative z-axis. 

The classical Hamiltonian per particle (see Appendix 
[a] ), expressed in terms of the canonical variables reads. 


hci{p,qJz^(^) 


{a,z\HD\a,z) 

j 



(4) 


The phase space of the classical Hamiltonian is IR^ x 
(for bosons and atoms respectively) and, except in some 
limiting cases (7 = 0 or cjo = 0)^ is non-integrable, having 
the energy per particle hcz(p, g, 0) = e as the only 
constant of motion. Both the QPT and the ESQPT are 
reflected in the classical energy surface. 

Since the number of bosons is not limited, the range 
of possible energies e is only lower bounded. The second- 
order QPT, according to the Ehrenfest classification, ap¬ 
pears as a discontinuity on the second derivative of the 
semiclassical ground state energy eo ( 7 ), which can be ex¬ 
pressed as mi ns US] 


eo(7) = 


-Wo for 7 < 7c, 

-‘f (^ + ^) for 7 > 7c- 


(5) 
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As the energy increases, in the superradiant region the 
energy surfaces acquire different structures, associated 
with the available phase space. They are marked by 
the ESQPT [131 ISH US]: sudden changes in the slope 
of the density of states. For energies in the inter¬ 
val e G [ 60 ( 7 ),—cJo] the surface of constant energy is 
formed by two disconnected lobes, which merges as the 
energy reaches e = —cjo- For energies in the interval 
e G [—the energy surface is formed by a sole 
lobe restricted to a fraction of the Bloch sphere. Finally, 
for energies larger than cjo? the whole Bloch sphere be¬ 
comes accessible. So, as mentioned, the QPT separates 
the system into a normal and a superradiant phase, with 
the latter having three regions separated by the ESQPTs 
iniiisiiiT], as well as a great richness of regularity and 
chaos. 


III. CLASSICAL AND QUANTUM DYNAMICS 
OF THE DICKE MODEL 


A. Poincare sections and Lyapunov exponents 


In order to determine the presence of regularity or 
chaos in the semi-classical system, we study the dynamics 
of the canonical variables. We employ Poincare sections 
for a qualitative insight and the Lyapunov exponent to 
quantify the presence of chaos [48l [49]. The surface of 
interest is defined by the plane in the variables 0 — jz 
which has p = 0 and satisfies Hci{q,p = O^jz^cj)) = e. 
This choice ensures a broad sampling of orbits because 
all of them intersects this surface. Under these conditions 
there are two different values of q±{jz, solutions 
of the quadratic equation Hci{q^p = 0 , 0 ) = e. 


27 / ~ 2 

(l±{izA,f) = - y^-jz COS(f) + 

UJ 


± 


1 47^ / _ 2\ 2 ~ 

y ^ q - j cos 2 0 + - (e - ^ojz). 


(6) 


The intersection of the classical orbits with this surface 
p = 0 defines the Poincare surface sections. 

Classical trajectories with energies very close to the 
ground state energy are regular. They can be de¬ 
scribed by an approximated quadratic, integrable Hamil¬ 
tonian, obtained by considering small oscillations around 
the minimum energy configuration [ 8 |. As the en¬ 
ergy increases the quadratic approximation breaks down. 
Chaotic trajectories appear when certain excitation en¬ 
ergy ech is reached, which is coupling dependent. For 
e > ec/i a region of soft chaos, characterized by a mix¬ 
ing of regular and chaotic orbits, is found. Fully clas¬ 
sical chaotic regions are always present at large enough 
excitation energies, both in the normal and superradi¬ 
ant phases m .These chaotic regions include the ground 
state in a small vicinity of 7 ^ 7 c, the QPT. A detailed 
study of regularity and chaos in the classical dynamics of 
the Dicke Hamiltonian, fully covering the regions of in¬ 



Figure 1: Poincare sections (top, (a) and (b)) and Lyapunov 
exponents (bottom, (c) and (d)) for the two sections q± as 
functions of jz and 0, for e = — 1.4ct;o and 7 = 27 c. In the 
Poincare sections the colors are associated with different clas¬ 
sical trajectories. For the Lyapunov exponents the color code 
is given on the bar. Blue depicts the regular regions. 


terest in energy and coupling constant will be presented 
elsewhere m- 

In Fig we present the Poincare sections and the Lya¬ 
punov exponents for e = — l.dc^o and 7 = 2jc in reso¬ 
nance {uj = cjo), for the two surfaces q±, as functions 
of jz and (j). Fig. [^a) and [^b) display the Poincare 
sections, with different colors for different orbits. For 
this energy, regular and chaotic regions coexist. Regu¬ 
lar regions are observed in Fig. [^a) at the top, in a 
semicircular area covering the whole available interval of 
variable 0. Another regular region can be identified at 
the triangular central sector, and a small third one at the 
bottom, around 0 = 0 and jz = —0.8. Similar regions 
are found in the other section in Fig. [^b). Figures [^c) 
and[^d) display the corresponding Lyapunov exponents, 
with their magnitudes represented by the color bar on 
the right. The blue regions are the regular ones (null 
Lyapunov exponent). It is apparent from the compari¬ 
son of the figures that the Lyapunov exponent quantifies 
the presence of chaos, which is qualitatively suggested by 
the Poincare sections. 


B. Quantum description of chaos 

In the previous studies [T 2 j [13], we have employed 
the Peres lattices m as quantum counterparts of the 
Poincare sections. They represent a qualitatively sen¬ 
sitive probe that allows to visualize the competition 
between regular and chaotic behavior in the quantum 
spectrum of a system [52]. We have also employed 
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the Anderson-Darling parameter to distinguish between 
Wigner and Poisson-like distributions of nearest neighbor 
energy differences. While useful, they can only be em¬ 
ployed to study energy intervals. In this work we move 
a step forward, employing two measures which allow to 
study the presence of regularity and chaos in the quan¬ 
tum regime, for every pointy in phase space at a given 
energy. They are the Husimi function and the Participa¬ 
tion Ratio of coherent states on the energy eigenstates. 


C. The Husimi function 

In order to both identify chaotic or regular character¬ 
istics in individual energy eigenstates, and to quantify 
chaos in the phase space, we use the Husimi function. 
The Husimi or Q function is one of the simplest distribu¬ 
tions of quasiprobability in phase space. Having a well- 
defined classical limit, it allows the comparison between 
the quantum and classical phase-space description of the 
Dicke model minis]- When j ^ oc the Husimi function 
reduces to a classical probability function on phase space 
obeying the Liouville equation mm- The Q function is 
defined as the expectation value of the density matrix in 
a set of coherent states. For the eigenstates \Ek) with 
energy they are the module squared of their projec¬ 
tion in the coherent state \a^z)^ given in Eq.(|^. The 
resulting function is 

Qk{a,z) = \{z,a\Ek)\^. (7) 

The Husimi function has been employed in the Dicke 
model by several authors to study the quantum-classical 
transition and equilibration mmis], the wave functions 
of individual states [43l [53] , and the ground-state QPT 

[HISS]. 

In order to make contact with the classical calculations, 
we evaluate the Q/c function along the same energy sur¬ 
faces, q± given in Eq.([^. Selecting a set of eigenstates 
whose eigenenergies satisfy \Ej^/j — e\ ~ 0. The com¬ 
parison between the different Husimi functions and the 
Poincare surface sections is quite illustrative. The evalu¬ 
ation of the Husimi function employing the efficient co¬ 
herent basis (see appendix]^ involves technical aspects, 
described in appendix [Cj 

Density plots of the Husimi function Q^{jz^<P) = 
QkUz,(l), g = g+,p = 0) for a j = 60 system with 7 = 2jc 
in resonance {uj = ujq) are shown in Eigj^ We chose five 
parity positive eigenstates, in the energy region of FigjlJ 
Elj = —IAuq. The similitude of the high density areas 
of the Husimi functions, depicted in white and light blue, 
with different Poincare sections at energy e = —IAujq, in¬ 
cluded at the center and superimposed in the quantum 
results, is noticeable. The shown eigenstates were chosen 
to reproduce the gross structures observed in the clas¬ 
sical results: the three largest stability islands (top and 
right panels) and the two largest chaotic seas in the left 
panels. 
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Figure 2: Density plots of the Husimi functions, Q^, for five 
parity positive Hamiltonian eigenstates of a j = 60 system 
with 7 = 27 c in resonance cjq = u. The energies of the 
five states are very close to E/j = — 1.4ct;o. The Husimi 
functions have strong resemblance to different Poincare sec¬ 
tions at energy e = — 1.4ct;o, which are exhibited at the 
center and superimposed in the quantum results. The in¬ 
dices and energies of the parity positive eigenstates (E^) 
are (from top in clockwise direction) = —1.4014ct;oj, 

E^gg = — 1.3967ct;oj, E^^ = — 1.3941ct;oj, E^g = —lA04:lujoj 
and E^gi = —lA064ujoj. 


IV. (DE)LOCALIZATION OF THE COHERENT 
STATES 

F. Haake suggest in na, in the context of the kicked 
top, that the minimum number of eigenstates of Flo- 
quet operators necessary to reconstruct a coherent state, 
Dmin^ can be a useful a tool to identify chaotic and reg¬ 
ular regimes. In [14] it is also shown that the scaling of 
Dmin with the dimension j (in our context to the num¬ 
ber of atoms), scales as ^/J in regular regions, and as j 
in chaotic ones. In the regular case the set of eigenstates 
tends to be localized, in correspondence with the classical 
regular movement inside stability islands. In the ther¬ 
modynamical limit j ^ 00 an infinitely small fraction of 
eigenstates (^ VI/j = l/VJ) is enough to reconstruct a 
coherent state associated with a regular region, whereas 
for a coherent state in a chaotic region this fraction goes 
to a finite value. This measure is proposed as an analogue 
of the classical Lyapunov exponent. 

That this localization in the space of eigenstates effec¬ 
tively takes place in the Dicke model can be seen in the 
distribution of the coherent state over the Hamiltonian 
eigenbasis = \ {a^ z\Ek)\‘^. These distributions 

are shown in Fig. |^and Fig. [^ for a regular and a chaotic 
point in phase space, respectively. Two different system 
sizes {j = 60 and 120) are shown. In both figures the dis¬ 
tributions become narrower as j increases, but the num¬ 
ber of participating states is clearly larger in the chaotic 
case respect to the regular one. In the bottom panels of 
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Figure 3: Top: Distribution of a coherent state over the en¬ 
ergy eigenstates, \C^(a,z)\‘^ = \{a, z\Ek)\‘^, for e = —IAujq, 
0 = 0 ? jz = —0.75, q = q+, p = 0 and j = 60,120. Bottom: 
3-D plot including {Jz)klj as a third coordinate for j = 60. 
The amplitudes are localized along a regular line. 


Figure 4: Same as Figbut for coherent parameters 0 = 0.0, 
jz = —0.55, p = 0 and ^ = ^+. The coherent state is mainly 
built upon eigenstates located in the disordered part of the 
Peres lattice (as seen in the bottom panel). 


the same figures, the distributions for the j = 60 cases 
are displayed in a 3-D plot against the energy eigenstates 
and {Jz)k/j- In the plane [E^jj^ {Jz)klj) a Peres lattice 
is formed. It has two regular arrays of points at the edges, 
associated with regular dynamics, and an interior region 
of scattered points which characterize the chaotic states 
m- In Fig. [^the intensities are clearly arranged along 
one regular edge, while in Fig. |^most of the intensity is 
located in the chaotic area. So, for the regular case, the 
coherent state is mainly built upon states in the ordered 
part of the Peres lattice, whereas for the chaotic case the 
main contribution comes from states in the disordered 
part. 

It should be emphasized that the numerical conver¬ 
gence is challenging when the numbers of atoms is in¬ 
creased. A careful analysis of this point is presented in 
appendix 


A. The participation ratio 

The use of the participation ratio Pr as a measure 
of localization of a quantum state was introduced several 
years ago [56ll6Q] . It has been applied to the Dicke model 
in connection to the study of equilibrium of many-body 
quantum closed systems m- Also, it has been employed 
to show that the equilibration process depends on the 
spreading of the initial state over the perturbed basis 
[62] . In this section we introduce the participation ra¬ 
tio Pr as a quantitative measure of the localization of 
the coherent states in the eigenstate basis. At variance 


with Dmin^ does not requires a cutoff (the smallest rel¬ 
evant contribution) to determine how many eigenstates 
are enough, the Pr has its own scale. 

For a pure quantum state |^), expanded in a basis 
{|0/e)} of dimension A, the participation ratio is, 


Pr 


1 


(8) 


It is defined in the interval Pr G When Pr = 1 

it means the state |^) is identical to one of the states 
of the basis, and it is considered as having maximum 
localization. On the other hand, if every state of the basis 
equally contribute to the state, we would have ( 0 n|^) = 
1/y/N. In this case, Pr = N. So, the maximum value 
of the Pr is related to maximum delocalization in the 
Hilbert space. 

In order to make contact with the classical phase space 
we employ the coherent states given in Eq.Q, whose 
parameters are defined by a single point {q^p^jz^(f) 
the phase space with energy e, as explained in Section 

III. 

The Pr is 


The Pr is obtained from the Husimi function, Eq. 
for every eigenstate evaluated over a single point in phase 
space. We restrict ourselves to cj = = 27c, and 

study the energy surface for a given e with p = 0. 

As a first example, we select the ground state en¬ 
ergy, and the point in phase space which characterize 
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Figure 5: 3-dimensional plots and respective density plots of 
Pr over the energy surface e = —IAujq for 7 = 27 c in res¬ 
onance uj = ujQ. Right and left columns correspond, respec¬ 
tively, to and q- with j = 60 and Nmax = 100 . 


the coherent state corresponding to the ground-state in 
the thermodynamic limit [21] [24] . Their canonical vari- 

ables are p = 0, g = ,and 

0 = 0 with energy as in Eq. (§. For 7 = 27c, it has 
Pr = 1.00585 for j = 30. It has nearly perfect overlap 
with the ground or the first excited state, which are de¬ 
generate in the superradiant phase (but still numerically 
distinguishable for j = 30). The value of Pr for this co¬ 
herent state becomes even closer to one for larger number 
of atoms. 

To explore the phase space for e = —IAujq, we calculate 
the Pr for points chosen over the energy surface using 
j = 60. The results are shown in Figure one for 
and the other for q_. We observe that the Pr over the 
surface closely resembles the distribution of Lyapunov 
exponents shown at the bottom of Fig. characterizing 
regions of regular or chaotic behavior. 

A detailed comparison between the Pr and the Lya¬ 
punov exponent is shown in Fig. where we plot points 
for q = q^ (top) and q = q- (bottom), using j = 60 
and Nmax = 100, along the line with 0 = 0. The Lya¬ 
punov exponent (black points) quantifies chaos, having 
non-zero values in chaotic regions and zero over regular 
ones. We can see there is a global agreement with the 
Pr (red points): lower values correspond to regular re¬ 
gions. However, this global agreement is clearer if we 
look into a binary criterion. We consider the quantity 
If Pi? < 1 we assign a zero value 
{PRbin = 0) just like in the Lyapunov exponent case (if 
A = 0, Ai)in = 0). For > 1 we assign the value one 
{PRbin = 1), as well as for a non zero Lyapunov exponent 
{Abin = !)• In Fig. 0 we show the results for the same 
482 points restricted to the 0 = 0 line. By only consid¬ 
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Figure 6 : Comparison between the Pr and Lyapunov expo¬ 
nent over the line e = —lAuo^p = 0,0 = 0. With j = 60 
and Nmax — 100. For q — qj^ (fop)? and q — q- (bottom). 
The black points stand for the Lyapunov exponent and the 
red ones for the Pr. 

ering the binary criterion the sensitivity of the Pr is re¬ 
markable. The global agreement for regular and chaotic 
regions is noticeable. The remaining differences can be 
attributed to the finite value of j. It can be expected 
that for larger j’s the regions with red and black points 
would perfectly match. 


B. Scaling of the Participation Ratio 

While we have shown that in the classical chaotic re¬ 
gions the Pr is larger than in the regular ones for given 
values of j and Nmax^ we still need to clarify if their 
magnitude is enough to determine if a point in the phase 
space is associated with regular or chaotic dynamics. In 
this section we demonstrate that its scaling as a function 
of j is a quantum measure of chaos associated with each 
point in phase space. 

Following the findings of Haake mi, we calculate the 
Pr as a function oiM = 2j for several points over chaotic 
regions or stability islands. The results are shown in 
Fig-i It can be seen that for points in regular regions, 
the Participation Ratio scales as VJV, while for points in 
chaotic regions it scales as Af. This confirms that the re¬ 
sults in [14] are also valid in the Dicke model, providing a 
second criterion for quantifying chaos employing purely 
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Figure 7: Same as Fig. [^but considering a binary criterion 
for the Pr and the Lyapunov exponent 

quantum tools. It follows that lim^^oo Pr/J^ goes to 
zero for a regular point, while for a chaotic point it re¬ 
mains constant, which corresponds to the binary criterion 
presented before. 


Figure 8: Black points: scaling of the Pr as a function of 
Af. The red lines are the fitted curved for each case, and the 
htting is shown in each figure. Calculated on a regular point 
determined by a zero Lyapunov exponent (top) and a chaotic 
point (bottom). Over the energy surface Eljuuo = —1.4 with 
p = 0 and 0 = 0. With a cut-off Nmax = 100 for every value 
of A/". In the insets we show the correspondent asymptotic 
constant values of Pr/VM and Pr/M for the regular and 
chaotic points, respectively. 


V. CONCLUSIONS 


regions of energies and coupling strengths is in progress. 


In the present work we have shown that it is possible 
to employ only quantum tools to characterize chaos in 
the phase space. We did that employing, in a qualitative 
way, the Husimi function to characterize the regular or 
chaotic behavior of individual eigenstates of the Dicke 
Hamiltonian. We showed that they can be associated 
with regular or chaotic sectors of the classical Poincare 
surface section. 

A quantitative quantum measure of chaos is provided 
by the Participation Ratio of the coherent state expanded 
in the basis of Hamiltonian eigenstates. In regular regions 
it scales as the square root of the number of atoms, and 
Pr/M smaller than 1, going to zero as M goes to infin¬ 
ity. In chaotic regions Pr/M tends to a constant, finite 
value. Its behavior was exhibited to follow closely that of 
the Lyapunov exponent. The Participation Ratio is, in 
this sense, the quantum equivalent of the Lyapunov ex¬ 
ponent, providing a quantum measure of chaos for each 
point of the phase space. A detailed study of the different 


Appendix A: The semiclassical Hamiltonian 


The state © expressed in the Fock and Dicke basis is, 
\a,z) = |a) 0 1^) = (Al) 

= E“.o Et.-, S I") ® li. ”•> ^ 


Then, the semiclassical Dicke Hamiltonian {az\H]j\az) 
(J = 1), reads 


/I - Izp 

Hci{a,z) =uj\a\'^ - ujQj P ^ 


+ 


7\/2j 


(A2) 


(a + (a*) {z -h 2 :*). 


1 + |.| 2 ' 

We introduced the real conjugate variables (g,p) trough 
a = JP ip) and {(p^jz) = (0, — cos6>) trough z = 























tan(6>/2)e*'^. By substituting the previous expressions in 
the expectation value of Hd!3 we obtain hd of Eq.Q. 


Appendix B: The efficient coherent basis 


There are two main challenges to solve the Dicke 
Hamiltonian. The Hamiltonian’s non-integrability im¬ 
plies it must be diagonalized numerically. Then, the 
second challenge is the dimension of the Hilbert space 
which is, formally, infinite. Therefore a cut-off procedure 
has to be implemented to numerically diagonalize and 
obtain the eigenvalues and eigenfunctions of the Dicke 
Hamiltonian. Instead of using the Fock basis to obtain 
the matrix elements, we employ the efficient coherent ba¬ 
sis (ECB), which is the exact Hamiltonian eigenbasis in 
the limit cjo = 0 0122]. The ECB is constructed from 
vacuum sates, |0)^^, of a new bosonic displaced operator 


Mt)JV 

(Bl) 

The displaced vacuum states are obtained from rotated 
(by — I around the ^-axis) atomic states, which means 
we are employing the eigenstates |j, of instead 
of those of Jz- The vacuum states of the new basis, ex¬ 
pressed in terms of the rotated raising SU(2) operator 

\^)mx ~ \^mx)\j^'^x) — 

\l U + 

where \o:mx) is a boson coherent state with Glauber pa¬ 
rameter arux = ~27ma,/(co’y27)- In previous works we 
have shown that employing this new basis with a given 
cut-off in the superradiant region, the number of con¬ 
verged states is orders of magnitude larger than those 
that can be obtained with the same cut-off in the Fock ba¬ 
sis. In other words, it is efficient in order to study larger 
systems or the ultra-strong coupling regime in compari¬ 
son with the Fock basis, where the diagonalization pro¬ 
cedure easily becomes intractable numerically. For more 
details about the efficient basis and convergence see 0 - 

Ej. 


Appendix C: The Husimi function in the efficient 
coherent basis 

We use the efficient basis (ECB) to diagonalize numer¬ 
ically the Dicke Hamiltonian. The k-th. eigenstate of the 
Dicke Hamiltonian |F/e), spanned in the ECB is (see Ap¬ 
pendix]^, reads 


\Ek) = {N;j,mx\Ek)\N-,j,mx), (Cl) 

N,mx 


where the coefficients = {N] j^mx\Ek) are calcu¬ 

lated numerically. The coherent state in Eq. (§ spanned 
in the eigenstate basis is. 


\a,z) = Y,{Ek\a,z)\Ek) = Y,C\a, z)\Ek). (C2) 

k k 

From the definition of the Husimi function, Eq. 0, it 
follows that Qk{oi^z) = \C^{a^ . The evaluation of 

the probability amplitudes 


C'"{a,z) = {Ek\a,z) = {N\j,mo,\a,z), 

N,m' 

(C3) 

involves the overlaps {N; j,mx\cx, z). By employing the 
definition of the ECB, we have 


{N-,j,m^\a,z) = a„^J^ \a,z) = 


^ {a-{amja){j,m::,\z), 


VM 


(C4) 


with = — 277713 ;The Glauber coherent 

states overlap is simply 

(a^Ja) = e-l“-^lV2e-l«lV2ea:„,a_ 


There is, however, an additional difficulty in estimating 
the overlap (j, mx\z)^ because \z) is a coherent state built 
on the basis of Jz, while the ECB is defined in terms of 
the eigenstates of |jf, m^,). To avoid the use of the 
Wigner D-matrix, we express the atomic coherent state 
in terms of the Jx eigenbasis. 


z) = Iw) = 


^ wJ - 


(l + kl^)' 


-j)x = 


2i 

( 1 + |M;| 2 )i 77 3 V 

Uix — J » 


E 


,J+^a 


b',’ 


(C6) 


From the previous expression it is direct to find the over¬ 
lap {jmx\z) = {3mx\w) in terms of the coherent param¬ 
eter w 


{jm:,\z) = ] \2\A ( ■ (C7) 

(1 + |m;|2)j y \j + rna;J 

To express the overlap in terms of the J;2-coherent param¬ 
eter z, we use the expectation values of the pseudospin 
operators. 


i^x) 

j 

{JA 

j 

{Jz) 


Re{z) \w\^ — 1 

\z\^ + 1 |u;p + 1 ’ 

Im{z) Im{w) 

\z\^ — 1 Re{w) 
\z\^ -h 1 |u;p + 1 ’ 


(C8) 


3 
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from which we obtain w{z) 

l-\- z 


(C9) 


The final expression for the coefficients C^{a,z)^ 
which give us the Husimi function Qf.(a, w(z)) = 

\C'‘{a,w{z))\'^, is 


C\a,z) = 

2^ \ \yN,mJ - ■ ■ - e e 2 X 


Nrrix 


2i 

j + rrix 


+ \wiz)yy 

27 ^ 

a-\ - ^m, 

U!y2j 


(CIO) 


C^’P{a,z) = {Ek\a,w{z)), which in terms of the BCE 
basis is, 


»(,))= E 

N,mx V + 


^rrix ,0) 

[{N; j, m^|a, w{z)) + p j, |a, w{z))] 


and whose overlaps have been already calculated. Fi¬ 
nally, after substituting the values of these overlaps, we 
obtain. 





\J-\-mx 


C^^P{a,w{z)) = 

X * , W\‘^Mo^rr 

^ 1 w{zy e 2 

) yjN\2{l+5E^ (1+I^(^)P)^' 

y (a - arux)^ 


X 


where w{z) is given by Eq. ( |C9[ ) and the coefficients 
C% are calculated numerically. 


Husimi function in the ECB with well defined parity. 


p { — 1)^ {aOLjnx)^ ^ 


With this, we have the Husimi function Q/c,p(<^,'?^(^)) = 
\C^^^{a,w{z))\^. 


We can take advantage of the Parity symmetry of the 
Dicke model to reduce, by a factor 2, the size of the 
Hamiltonian matrices that have to be diagonalized. The 
Hilbert space of the Dicke model has two invariant sub¬ 
spaces related each to the two eigenv alue s of the par¬ 
ity operator of Eq.([^. The ECB Eq.(C4) is not a set 
of eigenstates of the parity operator in Eq. Q, but it 
can be used to construct a basis with well defined parity 
(p = ±1). As it is shown in [12], the resulting basis is of 
the form, 


\N\j, mo,-,p) 


{\N]j,my) +p{-l)’^\N-,j, -my)) 


\/2(l + Sm^,o) 


(Cll) 


with A" = 0,1,... and > 0. 

Since the eigenfunctions of the Dicke Hamoltinian have 
well defined parity, they can expressed as follows, 

\Ek) = {N-,j,m’;p\Ek)\N-,j,m’-,p). (C12) 

N,m' 


Again, the = {N; j,m';p\E]^) are obtained numer¬ 

ically. In order to calculate the Husimi function we need 


Appendix D: Convergence of the numerical results 

We need to introduce a cutoff N^ax in the number of 
photonic-like excitations. The dimension of the space is 
Afst = For each coherent state, i.e. for 

each set of phase space parameters (a, z)^ we must guar¬ 
antee that it can be described in the truncated Hilbert 
space of the efficient basis, checking that its norm is close 
enough to one. The truncation limit the number Mconv 
of converged eigenstates, whose wave function is reliably 
described [Miaiisi. Additionally, we must check that 
the same coherent state can be described employing the 
converged eigenstates, satisfying 

conv 

\cHa,zy = i. (Di) 

k=l' 

As the number of atoms {2j) considered increase, the 
cut-off should be greater too. In Table |I| we list in the 
first four columns, the value of j, Nmax^-^st and Mconv 
It can be seen how fast the dimension of the space grows, 
and how the fraction of converged states diminish as j 
grows. 
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Norm. G. 
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only the converged states (with a bound in the accuracy of the wave function better than 10“^) for each value of j. Where 
A4t, J^conv, Norm., Norm. C. stand for ’Number of eigenstates’, ’Number of converged eigenstates’, ’Normalization’, and 
’Normalization considering only converged states’, respectively. 
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